home *** CD-ROM | disk | FTP | other *** search
/ CU Amiga Super CD-ROM 19 / CU Amiga Magazine's Super CD-ROM 19 (1998)(EMAP Images)(GB)[!][issue 1998-02].iso / CUCD / Programming / LEDA / source / src / graph_alg / _mc_matching.c < prev    next >
Encoding:
C/C++ Source or Header  |  1994-11-16  |  9.9 KB  |  464 lines
  1. /*******************************************************************************
  2. +
  3. +  LEDA  3.1c
  4. +
  5. +
  6. +  _mc_matching.c
  7. +
  8. +
  9. +  Copyright (c) 1994  by  Max-Planck-Institut fuer Informatik
  10. +  Im Stadtwald, 6600 Saarbruecken, FRG     
  11. +  All rights reserved.
  12. + 
  13. *******************************************************************************/
  14.  
  15.  
  16. /*******************************************************************************
  17. *
  18. *  MAX_CARD_MATCHING  (maximum cardinality matching)
  19. *
  20. *  An implementation of Edmonds' algorithm 
  21. *
  22. *
  23. *  J. Edmonds:  Paths, trees, and flowers 
  24. *               Canad. J. Math., Vol. 17, 1965, 449-467
  25. *
  26. *  R.E. Tarjan: Data Structures and Network Algorithms,
  27. *               CBMS-NFS Regional Conference Series in Applied Mathematics,
  28. *               Vol. 44, 1983
  29. *
  30. *
  31. *  running time: O(n*m*alpha(n))
  32. *
  33. *
  34. *  S. Naeher (October 1991)
  35. *
  36. *******************************************************************************/
  37.  
  38.  
  39. #include <LEDA/graph_alg.h>
  40. #include <LEDA/node_partition.h>
  41.  
  42.  
  43. enum LABEL {ODD, EVEN, UNREACHED};
  44.  
  45. static int cmp_degree(const node& v,const node& w)  // compare indeg of v and w
  46. { return indeg(v) - indeg(w); }
  47.  
  48.  
  49. static int greedy( graph &G, node_data<node>& mate )   // greedy heuristic
  50. { node v;
  51.   int count = 0;
  52.   forall_nodes(v,G)
  53.     if (mate[v] == nil)
  54.     { edge e;
  55.       forall_adj_edges(e,v)
  56.       { node w = target(e);
  57.         if (v != w && mate[w] == nil) 
  58.         { mate[v] = w;
  59.           mate[w] = v;
  60.           count++;
  61.           break;
  62.          }
  63.        }
  64.      }
  65.    return count;
  66.  }
  67.  
  68.  
  69.  
  70. static void heuristic( graph &G, node_data<node>& mate )
  71. {
  72.  // Heuristic (Markus Paul):
  73.  // finds almost all augmenting paths of length <=3 with two passes over
  74.  // the adjacency lists
  75.  // ("almost": discovery of a blossom {v,w,x,v} leads to a skip of the
  76.  // edge (x,v), even if the base v stays unmatched - it's not worth while
  77.  // to fix this problem)
  78.  // if all adjacent nodes w of v are matched, try to find an other 
  79.  // partner for mate[x], and match v and w on success
  80.  
  81.   node u,v;
  82.   bool found;
  83.   node_data<int> all_matched(G,false);
  84.  
  85.   forall_nodes( v, G )
  86.   {
  87.     edge e = G.first_adj_edge(v);
  88.  
  89.     if( mate[v]==nil )
  90.     {
  91.         while(e != nil && mate[target(e)]!=nil) e = G.adj_succ(e);
  92.  
  93.         if( e != nil ) 
  94.           { mate[v] = target(e);  
  95.             mate[target(e)] = v ;  
  96.            }
  97.         else // second pass
  98.           {
  99.         all_matched[v] = true;
  100.  
  101.             forall_adj_edges(e,v)
  102.             {
  103.               node w = target(e);
  104.               node x = mate[w];
  105.     
  106.               if( ! all_matched[x] ) 
  107.               {
  108.                 while((found=G.next_adj_node(u,x)) && (u==v || mate[u]!=nil));
  109.     
  110.                 if( found ) 
  111.                   { mate[u] = x;  mate[x] = u ;
  112.                     mate[v] = w;  mate[w] = v ;
  113.                     G.init_adj_iterator(x);
  114.                     break;
  115.                    }
  116.             else
  117.                    all_matched[x] = true;
  118.                }
  119.              }
  120.            }
  121.         } 
  122.     }
  123.  
  124. }
  125.  
  126.  
  127.  
  128.  
  129. static void find_path(list<node>& L, node_array<int>&  label,
  130.                                      node_array<node>& pred,
  131.                                      node_data<node>& mate,
  132.                                      node_array<edge>& bridge,
  133.                                      node x, node y)
  134. {
  135.   /* computes an even length alternating path from x to y begining with a 
  136.      matching edge (Tarjan: Data Structures and Network Algorithms, page 122)
  137.      Preconditions: 
  138.      a) x and y are even or shrinked
  139.      b) when x was made part of a blossom for the first time, y was a base 
  140.         and predecessor of the base of that blossom
  141.    */
  142.  
  143.  if (x==y) 
  144.  { // [ x ]
  145.    L.append(x);
  146.    return;
  147.   }
  148.  
  149.  if (label[x] == EVEN) 
  150.  { // [ x --> mate[x] --> path(pred[mate[x]],y) ]
  151.    find_path(L,label,pred,mate,bridge,pred[mate[x]],y);
  152.    L.push(mate[x]);
  153.    L.push(x);
  154.    return;
  155.  }
  156.  
  157.  if (label[x] == ODD) 
  158.  { // [ x --> REV(path(source(bridge),mate[x])) --> path(target(bridge),y)) ]
  159.    node v;
  160.    list<node> L1;
  161.    find_path(L,label,pred,mate,bridge,target(bridge[x]),y);
  162.    find_path(L1,label,pred,mate,bridge,source(bridge[x]),mate[x]);
  163.    forall(v,L1) L.push(v);
  164.    L.push(x);
  165.    return;
  166.   }
  167. }
  168.  
  169.  
  170.  
  171.  
  172. list<edge> MAX_CARD_MATCHING(graph& G, int heur)
  173. { 
  174.  
  175.  
  176.   // input: directed graph G, we make it undirected (symmetric) by inserting
  177.   //        reversal edges
  178.  
  179.  
  180.   list<edge> R;
  181.  
  182.   int    strue = 0;
  183.   bool   done  = false;
  184.  
  185.   node_array<int>  label(G,UNREACHED);
  186.   node_array<node> pred(G,nil);
  187.   node_array<int>  path1(G,0);
  188.   node_array<int>  path2(G,0);
  189.   node_array<edge> bridge(G,nil);
  190.  
  191.   node_data<node> mate(G,nil);
  192.  
  193.  
  194.   // make graph symmetric
  195.  
  196.   edge_array<edge> rev(G,nil);
  197.   compute_correspondence(G,rev);
  198.  
  199.   edge e = G.first_edge();
  200.  
  201.   while (e && rev[e]) e = G.succ_edge(e);
  202.  
  203.   edge last_edge = nil;
  204.   if (e !=nil)
  205.   { last_edge = G.new_edge(target(e),source(e));
  206.     R.append(last_edge);
  207.     while (e != last_edge)
  208.     { if (rev[e] == nil)
  209.       { edge e1 = G.new_edge(target(e),source(e));
  210.         rev[e] = e1;
  211.         R.append(e1);
  212.         }
  213.       e = G.succ_edge(e);
  214.      }
  215.    }
  216.  
  217.   edge_array<edge> reversal(G,nil);
  218.  
  219.   for(e = G.first_edge(); e != last_edge; e = G.succ_edge(e))
  220.   { edge r = rev[e];
  221.     reversal[e] = r;
  222.     reversal[r] = e;
  223.    }
  224.  
  225.  
  226. /*
  227.   edge_array<edge> reversal(G,nil);
  228.  
  229.   edge e = G.first_edge();
  230.   edge last_edge = G.last_edge();
  231.  
  232.   do { edge e1 = G.new_edge(target(e),source(e));
  233.        reversal[e] = e1;
  234.        reversal[e1] = e;
  235.        R.append(e1);
  236.        e = G.succ_edge(e);
  237.       }
  238.   while(e != last_edge);
  239. */
  240.  
  241.  
  242.  
  243.   switch (heur) {
  244.  
  245.   case 1: { greedy(G,mate);
  246.             break;
  247.           }
  248.  
  249.   case 2: { G.sort_nodes(cmp_degree);  // sort nodes by degree
  250.             heuristic(G,mate);
  251.             break;
  252.            }
  253.   }
  254.  
  255.  
  256.   while (! done)  /* main loop */
  257.   {
  258.     node_list Q;
  259.     node v;
  260.  
  261.     node_partition base(G);    // now base(v) = v for all nodes v
  262.  
  263.     done = true;
  264.    
  265.     forall_nodes(v,G)
  266.     { pred[v] = nil;
  267.   
  268.       if (mate[v] == nil)
  269.       { label[v] = EVEN; 
  270.         Q.append(v); 
  271.        }
  272.       else label[v] = UNREACHED;
  273.     }
  274.   
  275.     while (!Q.empty())    // search for augmenting path
  276.     {
  277.       node v = Q.pop();
  278.       edge e;
  279.  
  280.       forall_adj_edges(e,v)
  281.       {
  282.         node w = G.target(e);
  283.  
  284.         if (v == w) continue;    // ignore self-loops
  285.  
  286.         if (base(v) == base(w) || (label[w] == ODD && base(w) == w)) 
  287.         continue;   // do nothing
  288.  
  289.         if (label[w] == UNREACHED)
  290.         { label[w] = ODD;
  291.           pred[w] = v;
  292.           label[mate[w]] = EVEN;
  293.           Q.append(mate[w]);
  294.          }
  295.         else  // base(v) != base(w) && (label[w] == EVEN || base(w) !=w)
  296.         { 
  297.           node hv = base(v);
  298.           node hw = base(w);
  299.   
  300.           strue++;
  301.           path1[hv] = path2[hw] = strue;
  302.   
  303.           while ((path1[hw] != strue && path2[hv] != strue) &&
  304.                  (mate[hv] != nil || mate[hw] != nil) )
  305.           {
  306.             if (mate[hv] != nil)
  307.             { hv = base(pred[mate[hv]]);
  308.               path1[hv] = strue;
  309.              }
  310.   
  311.             if (mate[hw] != nil)
  312.             { hw = base(pred[mate[hw]]);
  313.               path2[hw] = strue;
  314.              }
  315.            }
  316.   
  317.           if (path1[hw] == strue || path2[hv] == strue)  // Shrink Blossom
  318.           { node x;
  319.             node y;
  320.             node b = (path1[hw] == strue) ? hw : hv;    // Base
  321.  
  322. #if defined(REPORT_BLOSSOMS)
  323.   cout << "SHRINK BLOSSOM\n";
  324.   cout << "bridge = "; 
  325.   G.print_edge(e);
  326.   newline;
  327.   cout << "base   = "; 
  328.   G.print_node(b);
  329.   newline;
  330.   cout << "path1  = ";
  331. #endif
  332.  
  333.             x = base(v);
  334.             while (x != b)
  335.             { 
  336.  
  337. #if defined(REPORT_BLOSSOMS)
  338.   G.print_node(x); 
  339.   cout << " ";
  340. #endif
  341.               base.union_blocks(x,b);
  342.               base.make_rep(b);
  343.  
  344.               x = mate[x];
  345.  
  346. #if defined(REPORT_BLOSSOMS)
  347.   G.print_node(x); 
  348.   cout << " ";
  349. #endif
  350.               y = base(pred[x]);
  351.               base.union_blocks(x,b);
  352.               base.make_rep(b);
  353.  
  354.               Q.append(x);
  355.  
  356.               bridge[x] = e;
  357.               x = y;
  358.              }
  359.  
  360. #if defined(REPORT_BLOSSOMS)
  361.   G.print_node(b);
  362.   newline;
  363.   cout << "path2  = ";
  364. #endif
  365.  
  366.             x = base(w);
  367.             while (x != b)
  368.             { 
  369.  
  370. #if defined(REPORT_BLOSSOMS)
  371.   G.print_node(x); 
  372.   cout << " ";
  373. #endif
  374.               base.union_blocks(x,b);
  375.               base.make_rep(b);
  376.  
  377.               x = mate[x];
  378.  
  379. #if defined(REPORT_BLOSSOMS)
  380.   G.print_node(x); 
  381.   cout << " ";
  382. #endif
  383.               y = base(pred[x]);
  384.  
  385.               base.union_blocks(x,b);
  386.               base.make_rep(b);
  387.  
  388.               Q.append(x);
  389.  
  390.               bridge[x] = reversal[e];
  391.  
  392.               x = y;
  393.              }
  394.  
  395. #if defined(REPORT_BLOSSOMS)
  396.   G.print_node(b);
  397.   newline;
  398.   newline;
  399. #endif
  400.           }
  401.           else  // augment path
  402.           {
  403.             list<node> P[2];
  404.  
  405.             find_path(P[0],label,pred,mate,bridge,v,hv);
  406.             P[0].push(w);
  407.             find_path(P[1],label,pred,mate,bridge,w,hw);
  408.             P[1].pop();
  409.  
  410.             for(int i=0; i<2; i++)
  411.               while(! P[i].empty())
  412.                { node a = P[i].pop();
  413.                  node b = P[i].pop();
  414.                  mate[a] = b;
  415.                  mate[b] = a;
  416.                 }
  417.  
  418.              Q.clear();                /* stop search (while Q not empty) */
  419.              done = false;             /* continue main loop              */
  420.              G.init_adj_iterator(v);
  421.              break;                    /* forall_adj_edges(e,v) ...       */
  422.            }
  423.  
  424.         } // base(v) != base(w) && (label[w] == EVEN || base(w) !=w)
  425.   
  426.       } // for all adjacent edges
  427.   
  428.     } // while Q not empty
  429.  
  430.   } // while not done
  431.  
  432.  
  433.  forall(e,R) G.del_edge(e);  // restore graph (only original edges in result!)
  434.  
  435.  list<edge> result;
  436.  
  437. /*
  438.   node v;
  439.   forall_nodes(v,G)
  440.     if (mate[v] != nil)
  441.     { edge e;
  442.       forall_adj_edges(e,v)
  443.         if (mate[target(e)] == v) 
  444.         { result.append(e);
  445.           break;
  446.          }
  447.      }
  448. */
  449.  
  450.  forall_edges(e,G) 
  451.  { node v = source(e);
  452.    node w = target(e);
  453.    if ( v != w  &&  mate[v] == w ) 
  454.    { result.append(e);
  455.      mate[v] = v;
  456.      mate[w] = w;
  457.     }
  458.   }
  459.  
  460.  return result;
  461.  
  462. }
  463.  
  464.